1. Persipan data dan library

2. Ambil data world.cities

data("world.cities", package = "maps")
set.seed(42)

cities <- world.cities %>%
  sample_n(300) %>%
  rename(lat = lat, long = long, pop = pop, name = name, country.etc = country.etc) %>%
  select(name, country.etc, lat, long, pop)

cities <- cities %>% filter(!is.na(pop) & pop > 0)

if(nrow(cities) < 200) {
  cities <- world.cities %>%
    filter(!is.na(pop) & pop > 0) %>%
    sample_n(300)
}

cities <- cities %>% mutate(logpop = log(pop + 1))
cities_scaled <- cities %>% select(lat, long, logpop) %>% scale()

knitr::kable(head(cities, 10))
name country.etc lat long pop logpop
Sarykamys Kazakhstan 47.78 78.72 4057 8.308445
Daosa India 26.88 76.33 76202 11.241156
Kaitangata New Zealand -46.25 169.85 860 6.758095
Cibeureum Indonesia -7.35 108.25 53443 10.886390
Suzano Brazil -23.55 -46.31 296003 12.598128
Deloraine Canada 49.18 -100.50 1078 6.983790
Santiago Mexico 25.41 -100.14 36498 10.505040
Hagavik Norway 60.18 5.43 1334 7.196687
Vardenis Armenia 40.18 45.72 11330 9.335298
Penarth UK 51.44 -3.17 23483 10.064075

Ringkasan dari 300 kota yang diambil dari data world.cities bawaan R

3. Menentukan jumlah cluster optimal

p1 <- fviz_nbclust(cities_scaled, kmeans, method = "wss") + labs(title = "Elbow Method - world.cities (sample)")

p2 <- fviz_nbclust(cities_scaled, kmeans, method = "silhouette") + labs(title = "Silhouette Method - world.cities (sample)")

subplot(ggplotly(p1), ggplotly(p2), nrows = 1, shareX = FALSE)
k_opt <- 3
cat('Pilihan k_opt =', k_opt, '\n')
## Pilihan k_opt = 3

Grafik menunjukkan metode Elbow maupun Silhouette sama-sama menyarankan k = 3 cluster optimal.

Metode clustering yang digunakan bersifat non-hierarkis (seperti K-Means), karena pembentukan kelompok dilakukan berdasarkan jumlah cluster (k) yang telah ditentukan sebelumnya, tanpa struktur bertingkat seperti pada metode hierarkis.

4. K-means clustering

set.seed(123)
kmeans_result <- kmeans(cities_scaled, centers = k_opt, nstart = 50)

cities_clustered <- cities %>% mutate(cluster = factor(kmeans_result$cluster))

cluster_sizes <- table(cities_clustered$cluster)
cluster_centers <- kmeans_result$centers

cat('Ukuran tiap cluster:\n')
## Ukuran tiap cluster:
print(cluster_sizes)
## 
##   1   2   3 
## 164  66  70
knitr::kable(head(select(cities_clustered, name, country.etc, pop, cluster), 12))
name country.etc pop cluster
Sarykamys Kazakhstan 4057 1
Daosa India 76202 2
Kaitangata New Zealand 860 2
Cibeureum Indonesia 53443 2
Suzano Brazil 296003 3
Deloraine Canada 1078 1
Santiago Mexico 36498 3
Hagavik Norway 1334 1
Vardenis Armenia 11330 1
Penarth UK 23483 1
Xacmaz Azerbaijan 37573 1
Vidisha India 143835 2

Tabel menunjukkan terbaginya kota menjadi 3 cluster berdasarkan populasi

Cluster 1: Kota kecil — populasi relatif rendah (ribu–belasan ribu).

Cluster 2: Kota menengah — puluhan hingga seratus ribu penduduk.

Cluster 3: Kota besar — populasi ratusan ribu ke atas.

164 kota di cluster 1

66 kota di cluster 2

70 kota di cluster 3

5. Visualisasi

Visualisasi clustering - Scatter (interaktif) + ggplot

p_scatter <- ggplot(cities_clustered, aes(x = long, y = lat, color = cluster, text = paste(name, country.etc, '<br>Pop:', pop))) +
  borders('world', colour = 'gray70', fill = 'gray95') +
  geom_point(aes(size = log(pop)), alpha = 0.8) +
  scale_color_brewer(palette = 'Set1') +
  labs(title = 'Peta Persebaran Kota Berdasarkan Cluster', x = 'Longitude', y = 'Latitude', color = 'Cluster', size = 'log(pop)') +
  coord_fixed(1.3) +
  theme_minimal()
## Warning: `borders()` was deprecated in ggplot2 4.0.0.
## ℹ Please use `annotation_borders()` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
ggplotly(p_scatter, tooltip = 'text')

Visualisasi menunjukkan bahwa hasil clustering ini memiliki validitas spasial — kota dalam cluster yang sama cenderung berdekatan di peta dunia.

Visualisasi interaktif dengan leaflet

pal <- colorFactor(palette = "Set1", domain = cities_clustered$cluster)

leaflet(cities_clustered) %>%
  addProviderTiles('CartoDB.Positron') %>%
  addCircleMarkers(~long, ~lat,
                   color = ~pal(cluster),
                   radius = ~scales::rescale(log(pop), to = c(3,10)),
                   stroke = FALSE, fillOpacity = 0.7,
                   label = ~paste0(name, ', ', country.etc, ' (Pop: ', format(pop, big.mark = ','), ')'),
                   popup = ~paste0('<b>', name, '</b><br>', country.etc, '<br>Pop: ', format(pop, big.mark = ','), '<br>Cluster: ', cluster)) %>%
  addLegend('bottomright', pal = pal, values = ~cluster, title = 'Cluster')

Visualisasi yang lebih interaktif

6. Evaluasi clustering

dist_euc <- dist(cities_scaled, method = "euclidean")

sil <- silhouette(kmeans_result$cluster, dist_euc)
sil_avg <- mean(sil[, 3])

cities_mat <- as.matrix(cities_scaled)

int_res <- intCriteria(cities_mat, as.integer(kmeans_result$cluster),
                       c("Calinski_Harabasz", "Dunn", "Davies_Bouldin"))

CH_index <- ifelse(!is.null(int_res$calinski_harabasz), int_res$calinski_harabasz, NA)
Dunn_index <- ifelse(!is.null(int_res$dunn), int_res$dunn, NA)
DB_index <- ifelse(!is.null(int_res$davies_bouldin), int_res$davies_bouldin, NA)

hasil_evaluasi <- data.frame(
  Silhouette_Avg = sil_avg,
  Calinski_Harabasz = CH_index,
  Dunn = Dunn_index,
  Davies_Bouldin = DB_index
)

hasil_evaluasi

Silhouette dan DBI menunjukkan clustering cukup baik, walau sebagian kota di perbatasan (populasi menengah) masih sedikit tumpang tindih antar cluster, ditunjukkan oleh nilai CH dan nilai Dunn.

7. Ringkasan dan penjelasan singkat

summary_tbl <- cities_clustered %>%
  group_by(cluster) %>%
  summarise(
    Rata_Populasi = round(mean(pop)),
    Max_Populasi = max(pop),
    Jumlah_Kota = n()
  ) %>%
  arrange(cluster)

knitr::kable(summary_tbl)
cluster Rata_Populasi Max_Populasi Jumlah_Kota
1 22179 292223 164
2 117411 1344436 66
3 76785 1123777 70

Interpretasi singkat:

  • Cluster 1: Kota kecil — populasi relatif rendah (ribu–belasan ribu).

  • Cluster 2: Kota menengah — puluhan hingga seratus ribu penduduk.

  • Cluster 3: Kota besar — populasi ratusan ribu ke atas.

Evaluasi Clustering:

  • Metode clustering bersifat non-hierarkis karena k sudah ditentukan di awal (k=3)

  • Nilai Silhouette rata-rata = 0.39 menunjukkan kualitas cluster cukup baik

  • Indeks Calinski–Harabasz = 149.6 menunjukkan konsistensi internal yang kuat, artinya variasi antar cluster lebih besar dibanding variasi dalam cluster.

  • Indeks Dunn = 0.0489 relatif kecil, menandakan jarak antar cluster masih bisa diperlebar.

  • Indeks Davies–Bouldin = 1.01 termasuk baik (mendekati 0), menunjukkan cluster yang cukup terpisah.

Kesimpulan: Struktur cluster sudah cukup optimal, namun masih dapat ditingkatkan dengan menyesuaikan jumlah cluster (k) atau metode normalisasi data.